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The paper focuses on the development of a general mathematical model 
and solution methodologies, to examine the behavior of thin structural 
elements such as beams, rings, and arches, subjected to large 
non- isothermal elasto-viscoplastic deformations. Thus, geometric as well 
as mater ial-type nonlinearities of higher order are present in the 
analysis. 

For this purpose a complete true abinito rate theory of kinematics and 
kinetics for thin bodies, without any restriction on the magnitude of the 
transformation is presented. A previously formulated elasto-thermo- 
viscoplastic material constitutive law is employed in the analysis. 

The methodology is demonstrated through three different straight and 
curved beams problems. Moreover importance of the inclusion of large 
strains is clearly demonstrated, through the chose applications. 
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Introduction 


The prediction of inelastic behavior of metallic materials at elevated 
temperatures has increased in importance in recent years. The operating 
conditions within the hot section of a rocket motor or a modern gas turbine 
engine present an extremely harsh thermo-mechanical environment. Large 
thermal transients are induced each time the engine is started or shut down. 
Additional thermal transients from an elevated ambient occur, whenever the 
engine power level is adjusted to meet flight requirements. The structural 
elements employed to construct such hot sections, as well as any engine 
components located therein, must be capable of withstanding such extreme 
conditions. Failure of a component would, due to the critical nature of 
the hot section, lead to an immediate and catastropic loss in power and 
thus cannot be tolerated. Consequently, assuring satisfactory long term 
performance for such components is a major concern for the designer. 

Traditionally, this requirement for long term durability has been a 
more significant concern for gas turbine engines rather than rocket motors. 
However, with the advent of reusable space vehicles, such as the Space 
Shuttle, the requirement to accurately predict future performance following 
repeated elevated temperature operations must now be extended to include 
the more extreme rocket motor application. 

Under this kind of severe loading conditions, the structural behavior 
is hightly nonlinear dut to the combined action of geometrical and physical 
nonlinearities. On one side, finite deformation in a stressed structure 
introduces nonlinear geometric effects. On the other side, physical 
nonlinearities arise even in small strain regimes, whereby inelastic 
phenomena play a particulary important role. From a theoretical 
standpoint, nonlinear constitutive equations should be applied only in 
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connection with nonlinear transformation measures (implying both 
deformation and rotations). However, in almost all of the works in this 
area (See Ref. 1), the two identified sources of nonlinearities are always 
separated. This separation yields, at one end of the spectrum, problems of 
large response, while at the other end, problems of viscous and/or 
non- isothermal behavior in the presence of small strain. 

The classical theories, in which the material response is 
characterized as a combination of distinct elastic, thermal, time 
independent inelastic (plastic) and time dependent inelastic (creep) 
deformation components cannot explain some phenomena, which can be observed 
in complex thermo-mechanical loading histories. This is particularly true 
when high-temperature non- isothermal processes should be taken into account. 
There is a sizeable body of literature^ »2 on phenomenological constitutive 
equations for the rate - and temperature - dependent plastic deformation 
behavior of metallic materials. However, almost all of these new "unified" 
theories are based on small strain theories and several suffer from some 
thermodynamic inconsistencies. 

In a previous paper3, the authors have presented an alternative 
constitutive law for elastic-thermo-viscoplastic behavior of metallic 
materials, in which the main features are: (a) unconstrained strain and 

deformation kinematics, (b) selection of reference space and configuration 
for the stress tensor, bearing in mind the rheologies of real 
materials, (c) an intrinsic relation which satisfies material objectivity, 
(d) thermodynamic consistency, and (e) proper choice of external and 
internal thermodynamic variables. Accuracy of the formulation was checked 
on a wide range of examples 1 *. 
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The formulation presented in this paper focuses on a mathematical 
model to examine the behavior of thin structural elements subjected to 
large non- isothermal elasto-viscoplastic deformations. Thus, geometric as 
well as material-type nonlinearities of higher order are present in the 
analyses. Such thin elements, including beams, rings and arches, are 
intended to present generic types of components, which might be located 
within or adjacent to the hot section of a rocket motor or gas turbine 
engine. 

The rate form of the constitutive equations suggests that^ a rate 
approach be taken toward the entire problem so that flow is viewed as 
history dependent process rather than an event. A direct consequence of 
the consistent adoption of the rate viewpoint in a spatial reference frame 
is that the problem is found to be governed by quasi-linear differential 
equations in time and in space. Hence, the analysis requires solution of 
an initial- and boundary- value problem involving instantaneously linear 
equations. The quasi-linear nature of the problem not only suggests an 
incremental approach to numerical solution, but also provides confidence in 
the completeness of the incremental equations. In this case, finite 
element solution capability is established; it should be noted, however, 
that the linearity of the instantaneous governing equations admits use of a 
wide variety of other established numerical procedures for spatial 
integration. A complete true ab initio rate theory of kinematics and 
kinetics for continuum and double curved thin structures, without any 
restriction on the magnitude of the strains or the deformation was 
formulated in Ref. 4 and will be rephrased here. 

Formulation of problems concerned with finite deformation of beams has 
followed two different paths**. Prescribing the beam by its deformed or 
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undeformed centroidal axis and cross section, one may introduce at the 
outset beam stress resultants and their conjugate kinematic variables 
characterizing displacement and rotation of the cross section. Together 
with appropriate beam constitutive equations and a global balance law a 
consistent theory is obtained. Alternatively, one may imbed beam theory in 
the setting of deformable solid continua, in which case one is concerned 
with local constitutive equations connecting the stress tensor with a 
strain tensor, which may in turn be expressed in terms of a combination of 
undetermined beam kinematic variables and functions of the beam coordinates. 
Momentum may then be balanced globally by integrating the local equations 
over the deformed beam configuration. Both paths will be considered in 
what follows. 

2. Two-Dimensional Plane Beams (A Plane Stress Problem) 

2.2 - Kinematics of the Continuum 

Let a continuum in space be described by two systems of coordinates, 

the x 1_ system, which stays at rest (the fixed system) and the u -system, 
which is associated with materials points (the convected material system). 

The transformation equations from one system to the other are: 

dx* = f * du° ( 1 ) 

a 

du a = u? dx* (2) 

l 

where 

(3 ) 

“ 3u° 


f a A 9u a 
1 = ax 1 


(4) 
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The covariant components of the metric tensor in the material system u° 
are: 


where Gjj 
system x*. 


g „ - f 1 G| 4 
00 o 0 ij 


(5) 


are the covariant components of the metric tensor in the fixed 
For the fixed cartesian system (Euclidian space) we have, 


Gij “ 6 ij (6) 

where are the components of the Kronecker delta. The coordinate lines 
of the x^-system are assumed to "deform" with the continuum in order to 

enable the material points to keep their coordinates (in the u°-system) 
unchanged. 

The contravar iant components of the velocity vector in the fixed 
system are defined by: 


v 


i 



( 7 ) 


It is impossible to define velocity as a change in the coordinates in the 
material system, however, distances are obviously changing. The length of 
the elementary arc in the material coordinates is given by: 

ds 2 = g du° du 8 (3) 

o0 

g 

Defining the rates of change in the material system by the rate of 

change of the elementary arc is, 


:(ds 2 ) 


H. 


at 1 '"*' 1 ~ 3t 
From Eq (9) one may conclude also that 


08 du a du 8 


i f a(ds) ^ 
ds ' dt ’ 


a f « i 9g o8 du a du 8 

M (log ds) • 2 — Si "as 


(9) 


( 10 ) 


The clue for the intrinsic rates of change may be unraveled, then, by the 
derivation of the rates of change of the metric tensor. 
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It can be shown 1 * that the rate of change of the metric tensor in the 
material coordinates is given by 


3g, 


a 8 

— — m V + v 
3t o ,8 8,0 


where 


v - g v Y = g f y v 1 » g f Y (f 1 ) 

a ,8 g Yo ,8 g Va i ,8 g Ya i 9t 1 6 ' 


(ID 


( 12 ) 


3 i i 

So — rr(fo) or v is the velocity in the fixed system as observed in 

ot p f p 

the material system. The components of the deformation rate tensor defined 
as follows, 


. A 1 9g a8 If i 

d o8 " 2 3t 2 V o, 8 4 V B,o^ ' d 8a 


and the components of the spin tensor as, 


03) 


^08 2^ V o,8 V 8,o^ 


8o 


Substitution of Eq. (13) into Eq. (10) yields, 


(14) 


-|j (log d3 ) - d oe ) 8 (15) 

As soon as the deformation rate is established as the time derivative 
of the metric tensor, the intrinsic characteristics of the continuum, being 
metric properties of space, are readily differented (with respect to time). 
For more details see Ref. 4. 

2.2 The Rate of Global Principles 

The principle of virtual power (or of virtual velocities), 


ij 


6 v 


J,1 


dV 


‘ { 
V 


p f J 6V. dV 


« 


vT J 6 v . dA = 0 


(16) 


is equivalent to the equations of equilibrium along with the complete set 
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of boundary conditions. In Eq. (16), are the contravariant components 
of the Cauchy stress tensor, p the mass density, and fJ is a vector of 
specific body forces. 

Total differentiation of Eq. (16) yields 1 *, 


1 l 2 3r*‘’ 1J <- v ‘- k " kJ)5v j.: 


dV 


f df J 

' j p “dt 


6v.dV 


f dT^ 

j Y *5t 6V j dA + 


r u fiii 

dt 


i dV - 


f J d6V J 

] pf ^ dV 




dA 


0 


(17) 


At any instant Eq. (17) must be satisfied. The virtual velocity and 
its time derivative are, then, independent. Moreover, the last three terms 
of Eq. (17) are equivalent to Eq. (16). Hence, the principle of the rate of 
virtual power may be obtained in its concise form. For further 
classifications, the total derivative of the stress components will be 
represented by the Jauman derivative, namely 



V c iJ 





( 18 ) 


and the following integrals are defined by 

V, 






! d - I <« 1J dk k - ° kJ d !k) 5 v j,i 

H 


dV 


V 

I 


ik r 

.k 0 6v ),i dv 


(19) 


(20) 


( 21 ) 
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Then, substitution in Eq. (17) yields the final form of the principle of 
the rate of virtual power, 

1 ■ 4 * l r " { " TE w j av * f 1 TE Sv j « (22) 

r A 

which is equivalent to 



i 



P 




= 0 


(23) 


and 


dlyT J -I 

dt 


d(a ij Vj) 

__ 


(24) 


A similar process can be applied 2 * to obtain the principle of the rate of 
balance of energy from the first law of thermodynamics. 

2.3 Constitutive Equations 

In a previous paper3, the authors have presented a complete set of 
constitutive relations for non isothermal, large strain, elasto-viscoplasti c 
behavior of metals. It was shown there3 that the metric tensor in the 
convected (material) coordinate system can be linearly decomposed into 
elastic and (visco) plastic parts. So a yield function was assumed, which 
is dependent on the rate of change of stress on the metric, on the 
temperature and a set of internal variables. Moreover, a hypoelastic law 
was chosen to describe the thermo-elastic part of the deformation. 

A time and temperature dependent "viscoplasticity" model was 
formulated in this convected material system to account for finite strains 
and rotations. The history and temperature dependence were incorporated 
through the introduction of internal variables. The choice of these 
variables, as well as their evolution, was motivated by thermodynamic 
considerations. 
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The nonisothermal elasto-viscoplastlc deformation process was 
described completely by "thermodynamic state" equations. Host 
investigators 1 * 2 (in the area of viscoplasticity) employ plastic strains as 
state variables. This study3 shows that, in general, use of plastic 
strains as state variables may lead to inconsistencies with regard to 
thermodynamic considerations. Furthermore, the approach and formulation 
employed in previous works leads to the condition that all the plastic work 
is completely dissipated. This, however, is in contradiction with 
experimental evidence, from which it emerges that part of the plastic work 
is used for producing residual stresses in the lattice, which, when 
phenomenologically considered, causes hardening. Both limitations were 
excluded from this formulation. 

The constitutive relation will be rephrased here as follows 


P 

a) if F = (tj " C p Q g p*)(tj - C p Q g gj) - k 2 (W.T) = 0 


(25) 


o*, and the temperature T 


where, s*, the Kirchhoff stress tensor, s* 
k ’ k 

are independent process variables, t* being the deviator of the Kirchhoff 


stress. s k> and W and B k are internal parameters, then 


with 


2Gv 


ui 


V 

1 + v 

1 




* « f s k * 

J 



c p„ g 




(26) 


l = 


11 


C p o 8 B k )(t l ~ C p o g fl2 ]V. . , ! 


(27) 
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( 28 ) 


b) 


and 


then 


with 


c) 




1 + 4 ni 


K - c 


P 0 8 4 1 


C P. 


S e 


i 
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P 

W 


- t i d k 

P o k i 


(29) 


V 

6 


i 

k 



(30) 


if 


F = 0 


(31 ) 



IE 

9T 


T > 0 


(32) 
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i 

k 



(33) 


P 

P V 

4 - 0 4' 2>( 4 - C|> o e s k> 


‘ ' d 2 |2(t ^ ' C P ° 8 % - TT f l <35) 


if F = 0 and -^y- s* + || T S 0 (36) 

9S k 

or F < 0 (37) 
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then 


(38) 



(39) 


(40) 


2. 4 Plane Stress Approximation 

By definition, a body is said to be in the state of plane stress 
parallel to the u 1 , u 2 plane when the stress components a 1 3, o 2 3, q 33 
vanish?. It is well known in literature that the case of plane stress is 
difficult to handle theoretically. Even linear elasticity has to treat 
this case in an approximate manner. To remove some of theoretical 
difficulties Durban and Baruch^ introduced the notion of Generalized Plane 
Stress, where instead of dealing with the quantities themselves, one deals 
with their average values. 

In our case the problem is even more difficult. The nonlinearities, 
which the general three-dimensional theory takes into account will also 
cause a large change of the geometrical quantities in the u3 direction. 
Clearly, some assumptions are needed to treat the case of plane stress as 
a two-dimensional case. 

The first basic assumption is that the thickness, h, of the plate 
defined by the coordinates u 1 , u 2 located in its middle plane, is small 
as compared with the other two dimensions. A second assumption is that the 
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external forces act in the u 1 , u 2 directions and are symmetrically 
distributed with respect to the middle plane. 

In a way similar to the procedure proposed by Durban and Baruch®, all 
the kinematic expressions are obtained by averaging the three-dimensional 
expressions. 

A basic assumption for the case of plane stress is that the components 
connected with the third direction are small and can be neglected. So, a 
new concept of generalized stress tensor is introduced 


i 

T k 


i u 

»k h 


(41) 


o 

It must be noted that in the linear theory of elasticity, where the 
geometry does not change, the averaged and generalized stress tensors 
coincide . 

So the three-dimensional incremental elasto-viscoplastic theory, 
developed previously, can be adopted for two-dimension plane stress 
problems. 


3. A Thin Curved Beam 

3.1 Doubly Curved Element 

A complete rate theory of kinematics and kinetics for doubly and 
singly curved thin structures, without any restriction on the magnitude of 
the strain or the deformation, was presented in Ref. 4. 

Five different shell theories (approximations), in rate form, starting 
with the simple Kirchhof f-Love theory and finishing with a completely 
unrestricted one, were considered there 1 *. 

The kinematic and kinetic equations for intrinsic shell dynamics, 
introduced in Ref. A are presented here, in compact form, together with 

basic notations. For simplicity we consider here Kirchhoff motion only®. 
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The components of the velocity vector w = y are 

o> a = y ' y o w° = y • y a , w n = y • n (42) 

where n is the unit normal to y, and y is the weighted motion. 

Expressions for the components of the velocity gradients, y, follow 
from differentiation of Eq. (42) 


oB 


o> „ - b „ w 
a; 8 aB n 


(43) 


w =» n 

a 


y „ m u + b . w 
8 n,a aB 


B 


(44) 


The time rates of the components of the metric and curvature tensors 
follow immediately from the above as 


d c8 + d Ba ’ ^aB w a;B d a ^B^ 


(45) 


To complete the kinematics, we get the components of the acceleration 
vector by time differentiation of Eq. (42) and through use of Eqs . (43) and 
(44), 


•a .a B ~o 
w ♦ d. w - o> w 
B n 


(46) 


y ♦ n = w + uw 
n 


The accelerations form the right sides of the equations of motion. The 
left sides are the static terms that can, for example, be expressed in 
terms of symmetrical stress resultants**. The result becomes 


•a — 1 r— 8 a ra 

W - rn [rr p + b A;g + m 


)S » 25“ . p“] - d° W B . 5“ w («7) 

A f p p n 
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w 
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- -1 r- aS r 
m n ' b 


aS 


-A — -a 6 

b b m 

3 aA 



p] 


-a 

u) 


w . 

a 


( 48 ) 


Here p 3 , p, m are loading components and mass, respectively, per unit area 

of y. 

3. > A Simplified Version of Curved Beam Element 



Fig. 1 - Reference Line of a Curved Beam 
A portion of the reference line for a curved beam is shown on Fig. 1. 
Toe current arc length is denoted by s, while <j» is the current angle of 
inclination of the normar to the reference line, and p is the radius of 

curvature. 

The stress resultants acting on the beam cross section are the bending 
moment M, the axial force N, and the shear force Q. The external load, 
measured per unit of current length of the reference line, has the components 

p and p in the direction of the unit vectors e and e respecti vely. 

3 n s n 
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If v g and denote the velocity components in the direction of the 


unit vectors e and e respectively, the rate of extension is 
s n 

d ■ 3s p 


(49) 


The rate of rotation, <t>, of a given section is given by 

3v v s 

w * * “ “aS + “ (50) 

while the generalized rate of deformation associated with bending is 

• 3* 3 , 3v n v s, 

k ' 15 - -35 ( -55 * ~> 151 > 

The rate of equilibrium equations for this simplified version may be 
put in the following form: 




|S , M »♦ . i |* , f , a P - o 
9s 3s 3s K n *n 


(52) 


3M 

3s 




Q + dQ 


0 


4. Numerical Solution 

The quasi-linear nature of the velocity equilibrium equations suggests 
the adoption of an incremental approach to numerical integration with 
respect to time. The availability of the field formulation provides 
assurance of the completeness of the incremental equations and allows the 
use of any convenient procedure for spatial integration over the domain B. 
In the present instance the choice has been made in favor of a simple first 
order expansion in time for the construction of incremental solutions from 
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the results of finite element spatial integration of the governing 
equations . 

The procedure employed permits the rates of the field formulation to 
be interpreted as increments in the numerical solution. This is 
particularly convenient for the construction of incremental boundary 
condition histories. 

The finite element method for spatial discretization has been well 
documented (see, e.g. Zienkiewicz^ or Oden 1 ^) and will not be detailed here. 
It should be noted, however, that as a consequence of the present 
formulation, the velocity equilibrium equations are not symmetric. This 
feature precludes implementation of a Ritz procedure as commonly employed 
in finite element analysis of infinitesimal deformation. Linear algebraic 
equations governing the discrete model for the finite case are developed 
employing the method of Galerkin. 

The spatial discretization results in: 

[K] iv} = { R} (53) 

where [K] is the nonsymmetric stiffness matrix, { v } is the vector 

containing the generalized nodal velocities, and {6} is the rate of the 
load. The solution to Eq. (53) at time to provides a basis for evaluation 

of a deformation increment and associated changes in internal stresses and 

boundary loading. The incremental solution defines the deformed 

configuration and stress rate at t = t Q + 5t thereby permitting definition 

of a new spatial problem at the later time. 

5. Applications 

The capabilities of the models presented here in have been evaluated 
through three simple numerical examples. The first example demonstrates 
the capability of the plane stress approximation to predict deflections and 
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ones in a beam loaded by a constant moment. Figure 2 illustrates the 
and tne °inite element model. A quarter of the beam was divided into 
elements in the vertical direction and into five elements in the 
octal nireetion. The external moment was introduced by six parallel 
s acting on tne section BC (see Fig. 2). 

Tne value of the external moment is 3500 kg/ cm, and the material of 
beam is C H <1 — 17. The viscoplastic properties of the material were 
neh experimentally from uniaxial tests in Ref. 1. This properties 
collaborated into the present material model. 

The variation of the deflection of point E as a function of time is 
i in Fig. 3. It is important to point out the value of the large 
mat ion analysis. After ten minutes of the deformation is increased by 
r.h at the same time there are important changes in the stress field 
"is* 4)* 



Fig. 2 - The Beam Model 
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The next example consists of a straight simply supported beam, loaded 
by a transverse concentrate force at the midspan. The beam is 25 inches 
long, two inches high and one inch wide. The material is stainless steel 
304 (Heat 9T2796). The material constants in sub section 2.3 were 
correlated with the uniaxial tension experimental results given in Ref. 12. 
The beam was subjected to a load of 2000 pounds at 1100°F, this load was 
then held constant for 312 hr., and then increased to 2250 pounds at 
1400°F. 

The primary purpose of this example is to compare the results, 
obtained by the two previously discussed models. The first one is the 
two-dimensional plane stress model, and the second one is the thin beam 
model as derived from thin shell theory. Figure 5 presents results in the 
form of load versus midspan deflection. The finite element model consists 
of five simple plane stress elements (dashed line in Fig. 5) or five 
sophisticated beam elements (full line in Fig. 5). 

It can be seen (Fig. 5) that the results agree quite well up to the 
312-hour hold period (points 3,4). During the hold period, the material 
hardens and only the beam model can represent this behavior after the load 
is further increased. 

The last example presents an analysis of a circular arch. The 
geometry of the shallow circular arch is shown on Fig. 6. The material is 
once again the 304 stainless steel. The arch is fixed at both ends and 
carries a concentrated load at the center. The elas to- v is co plas t i c 
analysis of this arch is performed with the aid of a ten curved beam elment 
model and with the inertia terms taken into account. The load P is assumed 
to be applied in a quasi-static manner at t = 0. The results of this 
analysis are shown on Fig. 6, as the time-history of the midspan 
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t w«r ('») 


Fig. 5 - Deflection vs Load 


ni<!sp.ir> ci ^placement , inches 



0 0.1 0.2 0.3 0.4 0.5 0.6 

time, hours 





Fig. 6 - Circular Shallow Arch 






displacement. The response of the arch starts with the instantaneous 
elastic deformation at t - 0, followed by slow deformation up to point B, 
whioh can be considered as a limit point for the given value of the load P. 


Beyond point B, 
may suggest the 


the displacements increase rapidly towards point C 
existence of critical time for the prescribed load. 


This 
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